#include <assert.h>
#include <stdio.h>
#include <stdlib.h>

#include "../src/distribution.c"
#include "../src/qfunc.c"

typedef struct _QFuncTestDatum QFuncTestDatum ;
struct _QFuncTestDatum 
{ 
  struct
    {
      float     kappa;
      float     logalpha;
      float     tau[(5-1)*4];
      float     loglambda[4];
      float     pz[5];
      Rating    ratings[8];
      float     logalphastddev;
    }                   inputs;
  struct
    {
      float q;
      float dqdlogalpha;
    }                   desired_output;
};

static void
test_qfunc (void)
{
  QFuncTestDatum data[] = 
    {{{0.511656,  0.939292, {0.453841, 0.24929, 0.474889, 0.8641, 0.794399, 0.415362,  0.569999, 0.920323, 0.588896, 0.1692, 0.127326, 0.947573,  0.177545, 0.607394, 0.276542, 0.661152}, {0.682189, 0.465102,  0.294869, 0.952559}, {0.880968, 0.891721, 0.783213, 0.013267,  0.427127}, {{3, 0}, {0, 4}, {1, 4}, {1, 0}, {0, 1}, {0, 4}, {2,  4}, {0, 0}}, 0.642432}, {-381.751, -371.923}}, {{0.308323,  0.149167, {0.632728, 0.227069, 0.738325, 0.228845, 0.0438323,  0.0578697, 0.610999, 0.281272, 0.866287, 0.450476, 0.334457,  0.62012, 0.184099, 0.985374, 0.0395883, 0.667561}, {0.303131,  0.0936529, 0.256375, 0.654294}, {0.876004, 0.451221, 0.948052,  0.505127,  0.243276}, {{2, 0}, {2, 2}, {2, 3}, {0, 4}, {3, 0}, {3, 0}, {1,  0}, {1, 1}}, 0.224152}, {-87.8321, -64.2376}}, {{0.209727,  0.276282, {0.199444, 0.166282, 0.598728, 0.99501, 0.333156,  0.715806, 0.264271, 0.37489, 0.149058, 0.730432, 0.224683,  0.707329, 0.845927, 0.636779, 0.968307, 0.0530351}, {0.969923,  0.185558, 0.0202553, 0.547908}, {0.726647, 0.961406, 0.810528,  0.271627,  0.527203}, {{0, 1}, {0, 4}, {2, 4}, {2, 4}, {2, 2}, {3, 2}, {0,  2}, {3, 4}}, 0.795124}, {-186.831, -166.982}}, {{0.2118,  0.276617, {0.194047, 0.0793173, 0.947529, 0.901727, 0.0449889,  0.348886, 0.722846, 0.194398, 0.199062, 0.712107, 0.754539,  0.141363, 0.229139, 0.526549, 0.734284, 0.593454}, {0.502492,  0.565143, 0.923756, 0.321827}, {0.975289, 0.77002, 0.711956,  0.0452106,  0.781242}, {{0, 3}, {2, 3}, {1, 4}, {0, 1}, {2, 2}, {3, 4}, {3,  2}, {2, 4}}, 0.690702}, {-216.152, -197.847}}, {{0.764427,  0.143484, {0.736254, 0.341817, 0.0415814, 0.949086, 0.537192,  0.62971, 0.287042, 0.807724, 0.308053, 0.103161, 0.552759,  0.21427, 0.805561, 0.538018, 0.629003, 0.892442}, {0.830272,  0.767998, 0.917047, 0.847232}, {0.0490291, 0.0772959, 0.15262,  0.703748,  0.312776}, {{3, 3}, {3, 4}, {0, 2}, {3, 3}, {3, 1}, {1, 2}, {0,  4}, {2, 2}}, 0.735479}, {-34.3038, -25.2102}}, {{0.111038,  0.754661, {0.775584, 0.105769, 0.823996, 0.946938, 0.467531,  0.00260762, 0.271237, 0.732668, 0.661971, 0.46459, 0.642234,  0.840226, 0.831699, 0.696591, 0.725187, 0.992994}, {0.78267,  0.619295, 0.572567, 0.289247}, {0.469894, 0.883816, 0.461528,  0.534585,  0.69431}, {{1, 3}, {2, 4}, {3, 2}, {1, 0}, {2, 3}, {1, 2}, {0,  1}, {3, 0}}, 0.778048}, {-185.036, -171.798}}, {{0.637532,  0.587648, {0.226779, 0.77544, 0.366295, 0.854979, 0.564808,  0.31085, 0.724062, 0.0147532, 0.733109, 0.614259, 0.998875,  0.0217587, 0.950439, 0.994963, 0.426309, 0.732512}, {0.480544,  0.111147, 0.96478, 0.197927}, {0.786234, 0.3331, 0.327248,  0.610279,  0.559455}, {{0, 1}, {1, 2}, {3, 1}, {2, 3}, {0, 3}, {3, 2}, {2,  3}, {1, 4}}, 0.55766}, {-149.566, -137.719}}, {{0.960952,  0.7553, {0.994647, 0.246809, 0.23689, 0.740547, 0.261537, 0.632551,  0.238015, 0.718788, 0.311098, 0.637587, 0.811707, 0.986276,  0.830554, 0.52644, 0.846926, 0.788349}, {0.0443201, 0.19334,  0.519679, 0.17807}, {0.484865, 0.635681, 0.558726, 0.42277,  0.490219}, {{3, 0}, {2, 1}, {0, 1}, {1, 1}, {2, 1}, {0, 4}, {0,  4}, {3, 4}}, 0.388871}, {-127.493, -113.044}}, {{0.321836,  0.682224, {0.228681, 0.756321, 0.0838208, 0.963436, 0.917583,  0.118734, 0.272114, 0.97716, 0.0870288, 0.592294, 0.425188,  0.188811, 0.0427087, 0.398953, 0.905509, 0.0107409}, {0.557843,  0.763273, 0.346782, 0.58797}, {0.0676248, 0.374401, 0.0249464,  0.905746,  0.838944}, {{3, 1}, {0, 1}, {3, 4}, {1, 4}, {0, 0}, {2, 2}, {2,  4}, {3, 0}}, 0.618081}, {-141.598, -131.846}}, {{0.941126,  0.94231, {0.921361, 0.499347, 0.669011, 0.96515, 0.834332,  0.907053, 0.243824, 0.776339, 0.791623, 0.5081, 0.338315,  0.765598, 0.23378, 0.744827, 0.991533, 0.177628}, {0.166155,  0.370426, 0.966586, 0.271881}, {0.327211, 0.752346, 0.0254606,  0.329571,  0.405851}, {{0, 3}, {0, 3}, {0, 4}, {2, 4}, {2, 2}, {1, 1}, {2,  2}, {2, 4}}, 0.252999}, {-222.46, -217.701}}} ;

  for (unsigned int i = 0; i < sizeof (data) / sizeof (data[0]); ++i)
    {
      ParameterizedNormalDistribution prior_logalpha = 
        parameterized_normal_distribution (data[i].inputs.logalphastddev);

      float hq = 0;

      for (unsigned int k = 0; k < 5; ++k)
        {
          hq -= data[i].inputs.pz[k] * log (data[i].inputs.pz[k]);
        }

      float dqdlogalpha;

      float q = qfunc (data[i].inputs.kappa,
                       data[i].inputs.logalpha,
                       &dqdlogalpha,
                       5,
                       data[i].inputs.tau,
                       data[i].inputs.loglambda,
                       data[i].inputs.pz,
                       data[i].inputs.ratings,
                       8,
                       &prior_logalpha.base,
                       hq);

      assert (fabs (q - data[i].desired_output.q) <=
              1e-4 * (1 + fabs (q) + fabs (data[i].desired_output.q)) ||
              (fprintf (stderr, "%g ?= %g\n", q, data[i].desired_output.q),
               0));

      assert (fabs (dqdlogalpha - data[i].desired_output.dqdlogalpha) <=
              5e-5 * (1 + fabs (dqdlogalpha) + fabs (data[i].desired_output.dqdlogalpha)) ||
              (fprintf (stderr, "%g ?= %g\n", dqdlogalpha, data[i].desired_output.dqdlogalpha),
               0));
    }
}

typedef struct _DQFuncTestDatum DQFuncTestDatum ;
struct _DQFuncTestDatum 
{ 
  struct
    {
      float     kappa;
      float     logalpha;
      float     tau[(5-1)*4];
      float     loglambda[4];
      float     pz[5];
      Rating    ratings[8];
      float     logalphastddev;
    }                   inputs;
  struct
    {
      float dtau[(5-1)*4];
      float dloglambda[4];
    }                   desired_output;
};

static void
test_dqfunc (void)
{
  DQFuncTestDatum data[] = 
    {{{0.356449,  0.364421, {0.571519, 0.345945, 0.112625, 0.588082, 0.779896,  0.837845, 0.77431, 0.822485, 0.546116, 0.093018, 0.782778,  0.644857, 0.379961, 0.722592, 0.816192, 0.372976}, {0.0527495,  0.970246, 0.790731, 0.043405}, {0.646899, 0.717248, 0.434282,  0.678984,  0.0753799}, {{1, 2}, {1, 4}, {3, 4}, {1, 0}, {3, 4}, {0, 2}, {0,  2}, {2, 4}},  0.371302}, {{-2.74099, -4.50504, 1.79125,  0.323143, -3.92378, -10.1928, -6.95732, -9.49896, -3.13673,  -4.50153, -7.23005, -8.00049, -2.55883, -4.99737, -6.64704,  -7.40505}, {-7.06356, -78.6456, -41.7014, -38.9482}}}, {{0.321657,  0.0909013, {0.295484, 0.533457, 0.547346, 0.268417, 0.749368,  0.440439, 0.764569, 0.62356, 0.369408, 0.717847, 0.948377,  0.250584, 0.316658, 0.747601, 0.157646, 0.207179}, {0.669759,  0.0303529, 0.723364, 0.528195}, {0.594379, 0.659051, 0.401707,  0.437293,  0.298895}, {{3, 3}, {2, 0}, {0, 1}, {1, 1}, {2, 0}, {3, 1}, {0,  0}, {1, 0}},  0.125594}, {{1.87011, 4.06616, 2.21765, 0.906883, 0.755418,  2.17733, 1.05318, 0.34759, 7.13659, 3.88955, 1.7753,  0.859951, -2.85277, -1.03873, -2.14226,  0.913173}, {-10.923, -4.48406, -14.3528, -12.6994}}}, {{0.854361,  0.168877, {0.549527, 0.685155, 0.089792, 0.545317, 0.180119,  0.967308, 0.141415, 0.294734, 0.863461, 0.219708, 0.983769,  0.0875549, 0.193702, 0.189355, 0.260406, 0.55936}, {0.599322,  0.530304, 0.858699, 0.122067}, {0.300427, 0.40471, 0.00433815,  0.95319,  0.7509}, {{0, 3}, {3, 2}, {2, 3}, {2, 0}, {2, 4}, {0, 0}, {0,  3}, {0, 2}},  0.719555}, {{1.08004, -1.14588, 1.62618, 3.54571, 0, 0, 0, 0,  1.87231, 0.768434, -5.68581, -1.91821, -0.454531, -0.905079,  1.6759, 0.584478}, {-32.5247,  0, -47.4215, -2.08019}}}, {{0.914546,  0.407873, {0.570781, 0.752247, 0.773131, 0.113139, 0.70732,  0.532539, 0.789362, 0.0255843, 0.513618, 0.343184, 0.528956,  0.466224, 0.914296, 0.812879, 0.670258, 0.344157}, {0.613869,  0.408169, 0.665919, 0.390968}, {0.862969, 0.688614, 0.751373,  0.983095,  0.292188}, {{2, 2}, {0, 3}, {0, 1}, {2, 0}, {3, 4}, {1, 3}, {0,  4}, {3, 2}},  0.936367}, {{-10.2509, -7.77146, -14.7678, -7.76828, -2.88856,  -4.51646, -6.58138, 0.665884, 3.60537, -0.25416, 4.58061,  1.02927, -6.34984, -9.68726, -4.95242, -7.01237}, {-76.0274,  -20.3159, -25.2895, -47.3026}}}, {{0.978242,  0.869956, {0.584868, 0.403829, 0.18888, 0.844371, 0.0712492,  0.0606449, 0.659924, 0.378147, 0.156953, 0.247766, 0.989667,  0.0339897, 0.543084, 0.839597, 0.323747, 0.643022}, {0.680115,  0.150983, 0.572374, 0.659927}, {0.387927, 0.214615, 0.594132,  0.789972,  0.803059}, {{2, 2}, {3, 2}, {2, 1}, {3, 1}, {1, 0}, {2, 2}, {1,  1}, {2, 3}},  0.810787}, {{0, 0, 0, 0, 5.70719, 11.5258, 6.16691,  2.64512, -6.85586, -0.908746, 3.2379, 11.9162, -4.69608, 2.68208,  12.4068, 2.54446}, {0, -37.3817, -49.1748, -27.4949}}}};

  for (unsigned int i = 0; i < sizeof (data) / sizeof (data[0]); ++i)
    {
      float dloglambda[4];
      float dtau[(5-1)*4];

      memset (dloglambda, 0, sizeof (dloglambda));
      memset (dtau, 0, sizeof (dtau));

      dqfunc (data[i].inputs.logalpha,
              5,
              data[i].inputs.tau,
              dtau,
              data[i].inputs.loglambda,
              dloglambda,
              data[i].inputs.pz,
              data[i].inputs.ratings,
              8);

      for (unsigned int k = 0; k < 4; ++k)
        {
          assert (fabs (dloglambda[k] - data[i].desired_output.dloglambda[k]) <=
                  1e-4 * (1 + fabs (dloglambda[k]) + fabs (data[i].desired_output.dloglambda[k])) ||
                  (fprintf (stderr, "%g ?= %g\n", dloglambda[k], data[i].desired_output.dloglambda[k]),
                   0));
        }

      for (unsigned int k = 0; k < (5-1)*4; ++k)
        {
          assert (fabs (dtau[k] - data[i].desired_output.dtau[k]) <=
                  1e-4 * (1 + fabs (dtau[k]) + fabs (data[i].desired_output.dtau[k])) ||
                  (fprintf (stderr, "%g ?= %g\n", dtau[k], data[i].desired_output.dtau[k]),
                   0));
        }
    }
}

int 
main (void) 
{
  srand48 (69);

  test_qfunc ();
  test_dqfunc ();

  return 0;
}
